Splitting of Long-period Comet C/2018 F4 (PANSTARRS)

Our images were photometrically calibrated using the methods described by Lister et al. (2022) with calviacat (Kelley & Lister 2019). Briefly, individual frames were calibrated to the Pan-STARRS 1 (PS1) photometric system (Tonry et al. 2012) using the ATLAS Refcat2 catalog (Tonry et al. 2018) and background stars. The photometric calibration includes color corrections to match the PS1 filters. Photometry of each nuclear component of C/2018 F4 was measured with circular apertures having a fixed linear radius of 8000 km projected at the distance of the comet (corresponding to angular diameters varying from ∼34 to 60, depending on the observer-centric distance of the comet). The background was determined within sufficiently large adjacent annuli, the size of which depended on the extent of the split comet in each image so as to get rid of contamination thereof as much as possible, but on the other hand, the annuli were set not to be too large so as to avoid contamination from background sources (annular radii range from 20''–120'').

We then proceeded to correct the photometric measurements for each nuclear component so as to remove the mutual contamination from the partially overlapping comae, simplistically assuming that each coma is in steady state (see Appendix for a detailed description of the derivation). In Figure 4, we plot the magnitude corrections for each nuclear component, both of which decrease with time as expected as the nuclear components apparently drifted farther apart from each other, resulting in diminishing mutual contamination. We show the apparent r-band light curves of the two components after the magnitude corrections in Figure 5, in comparison to the ones without correcting for mutual contamination. In general, the correction did not alter the overall trend of each component's apparent magnitude light curve, which overall became fainter over the course of the observed time span. At the beginning of our observing campaign, the apparent magnitude light curve of C/2018 F4-A was almost indistinguishable from that of C/2018 F4-B, but thereafter the latter appeared to drop slightly more rapidly, making it systematically fainter than the former (see Figure 5).

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

The apparent r-band magnitude of each nuclear component m_r is straightforwardly related to its total effective scattering cross section Ξ_e by

$$\begin{eqnarray}&&{{\Xi }}_{{\rm{e}}}=\displaystyle \frac{\pi }{{p}_{r}\phi \left(\alpha \right)}{\left(\displaystyle \frac{{r}_{{\rm{H}}}{\Delta }}{{r}_{\oplus }}\right)}^{2}{10}^{0.4\left({m}_{\odot ,r}-{m}_{r}\right)},\end{eqnarray} \tag{ 1 }$$

where r_H and Δ are the heliocentric and observer-centric distances, respectively, r_⊕ = 1 au is the mean Sun–Earth distance, p_r is the r-band geometric albedo, $\phi \left(\alpha \right)$ is the phase function normalized at zero phase angle, and m_⊙,r = − 26.93 is the apparent r-band magnitude of the Sun at r_⊕ = 1 au (Willmer 2018). We set p_r = 0.05 (Levasseur-Regourd et al. 2018) and adopted a linear phase model having a linear phase coefficient of β_α = 0.03 ± 0.01 mag deg⁻¹ (Meech & Jewitt 1987), both appropriate for comets. Here we did not attempt to incorporate the uncertainty in the geometric albedo, which cannot be constrained from our observations. Even if our assumed value is later found to be off, one can still easily scale and improve our estimates because of the inverse proportionality between the total effective scattering cross section and geometric albedo. We show the total effective scattering cross sections of both nuclear components varying with time in Figure 6, where the decline is obviously seen. We can think of two possibilities that may account for such trends—(1) cometary activity of both components faded as they receded from the Sun, and (2) smaller dust grains were more efficiently swept away by solar radiation pressure and drifted beyond the photometric apertures than were the larger counterparts. Based upon our analysis of the dust morphology (see Section 3.2), the former is preferred. Assuming a bulk mass density of ρ_d ∼ 1 g cm⁻³ for the dust grains, with the dominant dust grain size ${\bar{{\mathfrak{a}}}}_{{\rm{d}}}\sim 1\,\mathrm{mm}$ (see also Section 3.2), we estimated that the corresponding average net mass-loss rate over the course of our observing campaign was $\left\langle {\dot{M}}_{{\rm{d}}}\right\rangle =2{\rho }_{{\rm{d}}}{\bar{{\mathfrak{a}}}}_{{\rm{d}}}\left\langle {\dot{{\Xi }}}_{{\rm{e}}}\right\rangle /3\approx 9\pm 4\,{\rm{kg}}\,{{\rm{s}}}^{-1}$ for both nuclear components.

Zoom In Zoom Out Reset image size

Using our multiband observations, we examined the colors of the two nuclear components and the temporal trends thereof. Figure 7 is a comparison between the colors of Components A and B in terms of their g − r and r − i color indices. Unfortunately, our UH 2.2 m observation on 2020 November 20 was the only one in our campaign that covered the $i^{\prime}$ band, and therefore we could only compare the r − i color indices of the two components from the single night. What we found is that in the r − i regime, the mean color of Component A, ${\left\langle r-i\right\rangle }_{{\rm{A}}}=0.17\pm 0.06$, is not statistically different than the counterpart of Component B, ${\left\langle r-i\right\rangle }_{{\rm{B}}}=0.08\pm 0.06$, where the reported errors are both standard deviations from repeated measurements. As for the g − r regime, no statistically confident temporal color trend is seen for either of the nuclear components (see Figure 7). We then calculated the weighted mean values of their g − r colors over the course of our observing campaign to be ${\left\langle g-r\right\rangle }_{{\rm{A}}}=0.52\pm 0.04$ and ${\left\langle g-r\right\rangle }_{{\rm{B}}}=0.54\pm 0.04$ for Components A and B, respectively, which are again not statistically different. The errors thereof are standard deviations of the color measurements. Strictly speaking, what we measured is most likely the color from the dust environment rather than the nuclear components themselves, given the cometary appearances. Nevertheless, our measurements suggest that Components A and B had no major color differences.

Zoom In Zoom Out Reset image size

Licandro et al. (2019) measured the normalized spectral slope of C/2018 F4 before the splitting event to be ${S}^{{\prime} }=\left(4.0\pm 1.0\right) \%$ per 10³ Å from their spectroscopic observations. To compare pre- and post-split colors of the comet, we had to also calculate the spectral slope from the color measurements according to its definition by A'Hearn et al. (1984) and Jewitt & Meech (1986) as

$$\begin{eqnarray}&&{S}^{{\prime} }\left({\lambda }_{1},{\lambda }_{2}\right)=-\left(\displaystyle \frac{2}{{\rm{\Delta }}{\lambda }_{\mathrm{1,2}}}\right)\displaystyle \frac{{10}^{0.4\left[{\rm{\Delta }}{m}_{\mathrm{1,2}}-{\rm{\Delta }}{m}_{\mathrm{1,2}}^{\left(\odot \right)}\right]}-1}{{10}^{0.4\left[{\rm{\Delta }}{m}_{\mathrm{1,2}}-{\rm{\Delta }}{m}_{\mathrm{1,2}}^{\left(\odot \right)}\right]}+1},\end{eqnarray} \tag{ 2 }$$

where λ₁ and λ₂ mark the effective wavelength range, Δλ_1,2 ≡ λ₁ − λ₂ is the bandwidth, and Δm_1,2 and ${\rm{\Delta }}{m}_{1,2}^{\left(\odot \right)}$ are the color indices of the comet and the Sun, respectively. The normalization is chosen to be at the midpoint of the wavelength range. An object with ${S}^{{\prime} }=0$ corresponds to its color being completely the same as that of the Sun, while a negative value conveniently indicates a color bluer than the solar colors, otherwise redder. Substituting with the g − r and r − i color indices of the Sun and the corresponding bandwidths given in Willmer (2018), Equation (2) yields ${S}^{{\prime} }=\left(3.8\pm 3.0\right) \%$ per 10³ Å for Component A and $\left(4.8\pm 2.8\right) \%$ per 10³ Å for Component B as their normalized spectral slopes in the g − r regime, while their slopes in the r − i regime are $\left(3.3\pm 4.4\right) \%$ per 10³ Å and $\left(-2.5\pm 4.5\right) \%$ per 10³ Å, respectively. Thus, the colors of Components A and B were both consistent with that of pre-split C/2018 F4 reported by Licandro et al. (2019) within the noise level. This result, together with the similarity in the colors of the two components, may imply the homogeneity of the cometary nucleus of C/2018 F4. As a comparison, the two lobes of the nucleus of Jupiter-family comet 67P/Churyumov-Gerasimenko were also found to share a remarkably similar color (Sierks et al. 2015). The color of the comet we obtained resembles those of many other long-period comets reported in Jewitt (2015).

The dust morphology of comet C/2018 F4 implies the physical properties of its dust environment. Here, we applied the model by Finson & Probstein (1968) to study the dust morphology of the comet. Because the comet had split into two major fragments by the time of our earliest observation, each exhibiting a coma, we studied each component separately. In the Finson–Probstein model, all dust grains are assumed to leave the cometary nucleus at zero speed. The subsequent motion of the dust is then governed by parameter β, the ratio between accelerations of solar radiation pressure and solar gravitation, which is related to the grain size ${{\mathfrak{a}}}_{{\rm{d}}}$ and bulk density ρ_d by

$$\begin{eqnarray}&&\beta =\displaystyle \frac{{{ \mathcal C }}_{\mathrm{pr}}{{ \mathcal Q }}_{\mathrm{pr}}}{{\rho }_{{\rm{d}}}{{\mathfrak{a}}}_{{\rm{d}}}}.\end{eqnarray} \tag{ 3 }$$

Here, ${{ \mathcal C }}_{\mathrm{pr}}=1.19\times {10}^{-3}$ kg m⁻² is the proportionality constant, and ${{ \mathcal Q }}_{\mathrm{pr}}\approx 1$ is the scattering coefficient for typical cometary dust (Finson & Probstein 1968; Burns et al. 1979). Dust grains subjected to a common value of β but released over a range of different epochs form a syndyne, while those released at the same epoch from the cometary nucleus with different β form a synchrone.

We computed syndyne-synchrone diagrams for Components A and B of C/2018 F4, which were then compared with our observations (Figure 8). The orbital similarity of the two components renders no noticeable differences between their syndyne-synchrone diagrams. However, because of the close separation between the components and their relative positions, the dust feature of Component A around the position angle of the comet's negative heliocentric velocity projected onto the sky plane was strongly concealed by Component B, preventing us from meaningfully determining the dominant dust size of Component A; the exception to this is our last observation in 2021 December, when the two components were apparently much farther apart from each other but also much fainter than before. We also paid special attention to observations around the time when we were in the orbital plane of the comet and when we were the farthest from the plane. The former condition collapses all syndynes and synchrones into a single line and greatly helps us constrain the dust ejection speed, which is neglected in the Finson–Probstein model. The latter offers us a vantage point where different syndynes and synchrones would be maximally separated from our perspective.

Zoom In Zoom Out Reset image size

The general dust morphology of Components A and B remained largely unchanged in our observations, despite that the observing geometry, most notably the orbital plane angle, varied significantly (see Figure 1(d)). This indicates that the dust environment of the comet was dominated by large-sized particles. Indeed, comparisons between the syndyne-synchrone diagrams and observations consistently yielded that the main tail axis of Component B aligned nicely with syndynes of β ∼ 10⁻³ (corresponding to dust size ${{\mathfrak{a}}}_{{\rm{d}}}\sim 1$ mm) having release times no earlier than early 2020. As for Component A, using the last LOOK observation from 2021 December, we found that its dust morphology was also well matched with synchrones having β ∼ 10⁻³ and release epochs later than 2020 April, the same as what we obtained for Component B from the same observation. Judging from the nearly indistinguishable resemblance in the appearances of the two components around their respective inner regions throughout our observing campaign (see Figures 2 and 3), we infer that the dust environments of the two components were highly alike. Based on the dust release epochs, we deduced that the splitting of the cometary nucleus likely occurred in a time frame between early and mid-2020. Thereafter, both of the nuclear components were ejecting millimeter-sized dust grains in a protracted manner. The dominance of such large grains and the absence of small grains are qualitatively consistent with the model by Gundlach et al. (2015) where the small counterparts are retained by interparticle cohesion. Moreover, since the size of dominant dust grains remained unchanged, we favor that the observed decrease in the effective scattering cross sections of the two nuclear components (Figure 6) was more likely caused by their dwindling activity as the comet receded from the Sun, rather than due to dust grains being removed progressively from the photometric aperture by solar radiation pressure.

We now estimate the dust ejection speed, which was ignored in the syndyne-synchrone model. We noticed that isophotal contours in the inner regions of Components A and B remained largely circular throughout the observed period (see Figure 3, as an example), indicative of dominant dust grains being ejected largely isotropically at both nuclear components. The observation from UTC 2020 October 20 was taken closest to the plane-crossing time, granting us an unambiguous side view of the comet. From the observing geometry, we can derive that the out-of-plane width of a dust tail at such an epoch is related to the ejection speed of dust as

$$\begin{eqnarray}&&{V}_{{\rm{ej}}}=\displaystyle \frac{\tan {w}_{\perp }}{2{r}_{{\rm{H}}}}\sqrt{\displaystyle \frac{\beta {\mu }_{\odot }{\Delta }\sin \alpha }{2\tan \ell }}.\end{eqnarray} \tag{ 4 }$$

Here, V_ej is the ejection speed of dust, w_⊥ is the out-of-plane width of the dust tail at projected angular nucleus distance ℓ, and μ_⊙ = 3.96 × 10⁻¹⁴ au³ s⁻² is the heliocentric gravitational constant. From the observation, we measured the approximate width of the dust tail to be ${w}_{\perp }$ ≈ 35'' at nucleus distance ℓ ≈ 5'' for Component B. Substituting, we obtained V_ej ≈ 2 m s⁻¹. For Component A, the above method was inapplicable because of the contamination from Component B. Instead, we estimated the dust ejection speed of Component A using the turnaround distance of its dust coma in the sunward direction in the sky plane from

$$\begin{eqnarray}&&{V}_{{\rm{ej}}}=\displaystyle \frac{\sqrt{2\beta {\mu }_{\odot }{\Delta }\tan \ell \sin \alpha }}{{r}_{{\rm{H}}}},\end{eqnarray} \tag{ 5 }$$

where ℓ now refers to the angular nucleus distance of the turnaround point. After inserting the measured ℓ ≈ 20'' into the above equation, we find V_ej ≈ 2 m s⁻¹, which is similar to what we obtained for the dust ejection speed at Component B. Given the measurement uncertainties, these results are possibly no better than order-of-magnitude estimates. Nevertheless, the similarity in the ejection speeds of the dominant dust grains of the two nuclear components makes the two nuclear components appear to be even more homogeneous.

In order to explore the splitting event at comet C/2018 F4, we first measured the astrometry of both nuclear components with the Gaia Data Release 2 catalog (Gaia Collaboration 2018) as the reference for background sources. Although the nonsidereal tracking mode was turned on for our observations of the comet, field stars in the obtained images were not visibly elongated thanks to the apparent motion of the comet, image resolution, and seeing. Thus, we simply treated both nuclear components of the comet and background stars as bidimensional symmetric Gaussians to be fitted. The measurement uncertainties in the astrometric measurements of the nuclear components were properly propagated from errors in the astrometric reduction and centroiding. To extend the observed arcs of the components, we downloaded additional astrometry from the Minor Planet Center Database Search. ⁹

We then proceeded to compare the orbits of the two nuclear components by means of traditional orbit determination with FindOrb. ¹⁰ The astrometry obtained from the Minor Planet Center was debiased and weighted following the methods detailed in Eggl et al. (2020) and Vereš et al. (2017), respectively, because there were no available corresponding astrometric uncertainties. Our astrometry was simply weighted according to the measured astrometric uncertainties. Next, we employed FindOrb to fit the orbital elements of each nuclear component, taking into account perturbations from the major planets, Pluto, the Moon, and the 16 most massive asteroids in the main belt with the planetary and lunar ephemeris DE440 (Park et al. 2021), as well as relativistic effects. The observed-minus-calculated (O − C) astrometric residuals were largely consistent with the assigned or measured errors. Several measurements whose O − C residuals exceeded the 3σ level were downweighted accordingly. Nonetheless, this modification did not alter the initial orbital solutions beyond the 1σ uncertainty levels for both nuclear components. We list and compare the best-fitted orbital solutions in Table 1, where the orbital similarity is clearly manifested, strongly indicative of a common origin for the two components.

However, the orbital similarity of the two nuclear components itself is not sufficient to straightforwardly answer the questions of how and when the splitting event occurred. Therefore, we applied the fragmentation model first introduced by Sekanina (1977, 1978), which simplifies a splitting event as an instantaneous separation between two nuclear components at some fragmentation epoch, with a component leaving the other one with some relative velocity. The model neglects mutual gravitational interactions between the two components. Despite the similarity in the orbits of the two nuclear components, their orbital elements are statistically different, often indicative of a nonzero relative speed between the components at the time of splitting.

We employed our code, which incorporates the Levenberg–Marquardt optimization code MPFIT (Markwardt 2009), and the same gravitational perturbers with DE440 alongside relativistic corrections were included, to study the splitting at comet C/2018 F4. Previously we had exploited the code to delve into the fragmentation events at active asteroids P/2016 J1 (PANSTARRS) and 331P/Gibbs (Hui et al. 2017; Hui & Jewitt 2022). We first regarded C/2018 F4-A as the primary component of the comet and fitted the relative astrometry of C/2018 F4-B with respect to the former. In cases of no simultaneous measurements and only the astrometry of Component B being available, the relative astrometry was computed together with the obtained orbital elements of C/2018 F4-A listed in Table 1. Initially, we attempted to find a best-fit solution to the fragmentation epoch t_frg and separation velocity decomposed into radial, transverse, and normal (RTN) components referenced to the primary that would minimize the goodness of fit for the whole observed arc of Component B spanning from 2020 September 12 to 2021 December 2. However, regardless of our initial guessed values, the solution (here termed Solution I of the gravity-only model) always converged to a statistically similar result containing an unacceptably strong systematic trend in the O − C astrometric residuals for observations starting from 2021 August when the comet became observable again after the solar conjunction in 2021 (see Figure 9(a)). We then fitted the observations before and after the solar conjunction separately, thereby obtaining much better solutions (termed Solutions II and III, respectively, of the gravity-only model) without any conspicuous systematic trend in the O − C astrometric residuals, characterized by their goodness-of-fit values both being clearly smaller than that of Solution I (see Table 2). While separating the astrometry into two distinct arcs helped circumvent the systematic trend in the astrometric residuals, we are aware that the normal components of the best-fitted separation velocities in Solutions II and III do not agree statistically with each other, suggestive of problems in the gravity-only fragmentation model.

Zoom In Zoom Out Reset image size

Table 2. Best-fit Fragmentation Models of C/2018 F4 (PANSTARRS): Separation of Component B from Component A

Quantity		Gravity-only Model			Nongravitational Model
		Solution I	Solution II	Solution III	$\sim {r}_{{\rm{H}}}^{-2}$
Fragmentation epoch (TDB) a	tfrg	2020 May 2.8 ± 1.3	2020 Apr 26.8 ± 1.7	2020 Jan 23 ± 84	2020 Apr 22.3 ± 8.3
Separation velocity (m s−1)
	ΔVR	+3.145 ± 0.015	+3.031 ± 0.022	+2.45 ± 0.75	+2.93 ± 0.16
	ΔVT	−0.577 ± 0.023	−0.500 ± 0.025	+0.41 ± 0.62	−0.52 ± 0.17
	ΔVN	−0.1526 ± 0.0023	−0.1562 ± 0.0024	+0.143 ± 0.015	−0.363 ± 0.011
Differential NG parameter (au d−2) b
	ΔAR	0	0	0	$\left(+2.69\pm 0.37\right)\times {10}^{-8}$
	ΔAT	0	0	0	$\left(+1.17\pm 0.73\right)\times {10}^{-8}$
	ΔAN	0	0	0	$\left(+2.04\pm 0.13\right)\times {10}^{-8}$
Observed arc		2020 Sep 12-2021 Dec 02	2020 Sep 12-2021 Feb 06	2021 Aug 14-2021 Dec 02	2020 Sep 12-2021 Dec 02
Number of observations used		192	173	19	192
Normalized residual rms		1.364	0.804	0.755	0.710

Notes. The reported uncertainties are 1σ formal errors. Here, we treat Component A as the primary nucleus and fitted the astrometric observations of Component B. In each solution, the obtained separation velocity and differential nongravitational parameter of Component B decomposed into the RTN directions are referenced with respect to Component A at the corresponding fragmentation epoch. In the gravity-only model, the differential nongravitational parameter is held fixed at zero.

^a The corresponding uncertainties are in days. ^b Differential nongravitational parameter.

Download table as:  ASCIITypeset image

Thus, we decided to further take into account nongravitational effects, which may arise from anisotropic mass-loss activity. Because the fragmentation model in essence requires relative positions between the secondary and the primary components, technically the nongravitational acceleration yielded by the model is referenced with respect to the primary. Given the orbit of the comet (perihelion distance q = 3.4 au), we do not posit that nongravitational effects at C/2018 F4 would be driven by sublimation of water ice, because of the low equilibrium temperature. Therefore, the nongravitational force model by Marsden et al. (1973) based on isothermal water-ice sublimation is most likely inapplicable. Rather, we adopted a nongravitational force model varying with ${r}_{{\rm{H}}}^{-2}$, following Sekanina (1977, 1978), who incorporated nongravitational effects in the radial direction only. However, we found through testing that adding a radial differential nongravitational acceleration to be fitted was insufficient to completely remove the strong systematic trend in the O − C astrometric residuals for the full observed arc, although the goodness of fit was improved. Thus, we further incorporated the transverse and normal components of the differential nongravitational acceleration as additional free parameters in the fragmentation model, whereby the systematic trend was eliminated successfully, with the O − C residuals all consistent with the measurement uncertainties (see Figure 9(a)). We tabulate the best-fit parameters of the nongravitational fragmentation model in Table 2.

To demonstrate the reliability and robustness of the results, we repeated our aforementioned steps but instead treated C/2018 F4-B as the primary and fitted the relative astrometry of C/2018 F4-A. We again had to fit the observations separately with the gravity-only model; otherwise, the resulting O − C astrometric residuals in the post-conjunction observations would contain a strong systematic trend well beyond the noise level (see Figure 9(b)). Thus, the nongravitational model was invoked, with which we obtained almost identical best-fit results, except that the signs of RTN components of the separation velocity and differential nongravitational parameter are flipped, because here the roles between the two nuclear components are swapped (see Table 3).

Table 3. Best-fit Fragmentation Models of C/2018 F4 (PANSTARRS): Separation of Component A from Component B

Quantity		Gravity-only Model			Nongravitational Model
		Solution I	Solution II	Solution III	$\sim {r}_{{\rm{H}}}^{-2}$
Fragmentation epoch (TDB)	tfrg	2020 May 2.8 ± 1.3	2020 Apr 26.6 ± 1.6	2020 Jan 23 ± 82	2020 Apr 22.5 ± 8.1
Separation velocity (m s−1)
	ΔVR	−3.146 ± 0.015	−3.028 ± 0.022	−2.45 ± 0.74	−2.94 ± 0.15
	ΔVT	+0.575 ± 0.022	+0.497 ± 0.025	−0.43 ± 0.60	+0.53 ± 0.16
	ΔVN	+0.1522 ± 0.0023	+0.1558 ± 0.0023	−0.143 ± 0.015	+0.363 ± 0.011
Differential NG parameter (au d−2)
	ΔAR	0	0	0	$\left(-2.72\pm 0.36\right)\times {10}^{-8}$
	ΔAT	0	0	0	$\left(-1.20\pm 0.72\right)\times {10}^{-8}$
	ΔAN	0	0	0	$\left(-2.04\pm 0.13\right)\times {10}^{-8}$
Observed arc		2020 Sep 12-2021 Dec 02	2020 Sep 12-2021 Feb 06	2021 Aug 14-2021 Dec 02	2020 Sep 12-2021 Dec 02
Number of observations used		199	178	21	199
Normalized residual rms		1.374	0.834	0.722	0.734

Note. Same as Table 2, except that Component B is now treated instead as the primary, and we fitted the astrometry of Component A. In each solution, the separation velocity and the differential nongravitational parameter of Component A are expressed in terms of the RTN components relative to Component B at the instant of splitting.

Download table as:  ASCIITypeset image

In addition, we were concerned about how tailward biases in our astrometric measurements of the two components of comet C/2018 F4, if present, would influence the best-fit results. Such biases have been previously observed to exist for astrometric measurements of comets whose comae are asymmetric due to solar radiation pressure and/or anisotropic cometary activity, and therefore extrapolating to zero-aperture astrometry would be desirable (e.g., Hui & Ye 2020). However, because neither of the nuclear components of C/2018 F4 had sufficiently high signal-to-noise ratios, nor appeared strongly elongated (see Figures 2 and 3), we found that the zero-aperture astrometry had enormous uncertainties often exceeding ∼1''. Therefore, we saw no benefit in using zero-aperture astrometry. Rather, we utilized astrometry obtained with different aperture sizes, only to find that the nongravitational model would be needed to fully eliminate the systematic trend in the O − C astrometric residuals, and that the obtained best-fit parameters are all consistent within the noise level. Now we can firmly conclude that comet C/2018 F4 split into two major components most likely in late 2020 April, broadly consistent with the result from our syndyne-synchrone computation (Section 3.2), at a heliocentric distance of ∼3.7 au with a separation speed of 3.00 ± 0.18 m s⁻¹ between the two components, mostly in the radial direction. Furthermore, the best-fit models reveal that Component B was subject to a nongravitational acceleration stronger than that of Component A, statistically important in a plane perpendicular to the orbital transverse direction. Here, we feel the necessity to comment on the detection of the differential nongravitational effect between the two nuclear components of the comet.

At first glance, the detection of the differential nongravitational acceleration may seem to contradict our results from traditional orbit determination where gravity-only solutions were found to be sufficient, but in actuality, it does not. The reason is that the common origin of Components A and B is completely neglected in the traditional orbit determination, and yet this turns out to be the most important constraint of all. If we also incorporate the nongravitational model in FindOrb, there will be no statistically confident detection of the nongravitational parameters except in the normal direction for C/2018 F4-B, ${A}_{{\rm{N}}}=\left(+2.01\pm 0.50\right)\times {10}^{-8}$ au d⁻². However, the resulting differential nongravitational parameter in the normal direction is still statistically insignificant, as the uncertainty of the other component is dominant. Nonetheless, the result is still in agreement with that from the nongravitational fragmentation model within the 1σ level. We thus believe that the differential nongravitational acceleration between the two components of C/2018 F4 is authentic. This finding possibly indicates a smaller size of Component B compared to Component A of C/2018 F4, but it is also possible that this was a consequence of the mass loss at Component B being more anisotropic.